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We examine the thermal evolution of a sequence of compact objects containing low-mass hadronic 
and high-mass quark-hadronic stars constructed from a microscopically motivated equation of state. 
The dependence of the cooling tracks in the temperature versus age plane is studied on the variations 
of the gaplessness parameter (the ratio of the pairing gap for red-green quarks to the electron 
chemical potential) and the magnitude of blue quark gap. The pairing in the red-green channel 
is modeled assuming an inhomogeneous superconducting phase to avoid tachionic instabilities and 
anomalies in the specific heat; the blue colored condensate is modeled as a Bardeen-Cooper-Schrieffer 
(BCS)-type color superconductor. We find that massive stars containing quark matter cool faster in 
the neutrino-cooling era if one of the colors (blue) is unpaired and/or the remaining colors (red-green) 
are paired in a inhomogeneous gapless superconducting state. The cooling curves show significant 
variations along the sequence, as the mass (or the central density) of the models is varied. This 
feature provides a handle for fine-tuning the models to fit the data on the surface temperatures of 
same-age neutron stars. In the late-time photon cooling era we observe inversion in the temperature 
arrangement of models, i.e., stars experiencing fast neutrino cooling are asymptotically hotter than 
their slowly cooling counterparts. 



I. INTRODUCTION 

After the initial phase of rapid (of order of days to 
weeks) cooling from temperatures T ~ 50 MeV down to 
0.1 MeV, a neutron star settles in a thermal quasiequi- 
librium state which evolves slowly over the time scales 
10^ - 10^ yr down to temperatures T - 0.01 MeV. The 
cooling rate of the star during this period is determined 
by the processes of neutrino emission from dense mat- 
ter, whereby the neutrinos, once produced, leave the star 
without further interactions. The understanding of the 
cooling processes that take place during this neutrino ra- 
diation era is crucial for the interpretation of the data 
on surface temperatures of neutron stars. While the long 
term features of the thermal evolution of neutron stars 
are insensitive to the initial rapid cooling stage {t < 
years), the subsequent route in the temperature versus 
time plane, which includes in addition to the neutrino 
emission era the late-time {t > 10^ yr) photoemission 
era, strongly depends on the emissivity of matter during 
the neutrino-cooling era [ll-|3| • 

The neutrino emission during the neutrino-cooling era 
chiefly originates from the core of a compact object, i.e., 
from matter at and above densities ~ 0.16 fm~^ and is 
sensitive to its composition. Hadronic matter emits neu- 
trinos mainly via the beta and inverse beta decay pro- 
cesses and the charge neutral processes of Cooper pair 
breaking. The former are suppressed by pairing corre- 
lations {i.e., the emergence of a gap in the quasiparticle 
spectrum separating the ground state from the excited 
states below certain critical temperature), whereas the 
latter are the consequence of the formation of Cooper 
pair condensates [![. 

The densities in the centers of massive compact ob- 
jects can be by an order magnitude greater than the 
nuclear saturation density. The quark substructure of 
the baryons will play an increasingly important role as 



larger densities are reached in progressively massive ob- 
jects. The recent timing observations of binary millisec- 
ond pulsar J1614-2230 imply that compact objects with 
masses of order of 2 solar masses exist in nature 

It is the purpose of this work to examine the cool- 
ing behavior of massive compact objects containing cores 
of deconfined interacting quark matter. Specifically, we 
shall assume that above certain density the quarks are in 
continuum states filling a Fermi sphere. The mechanism 
by which the neutrinos are produced in quark matter are 
basically the same as in the baryonic matter, i.e., beta 
decays of d quarks and inverse beta decays via electron 
capture on a u quark. Because of our assumption of filled 
Fermi spheres and attractive interaction among quarks 
in dominant color antisymmetric states, the pairing of 
quarks - known as color superconductivity - needs to be 
taken into account In this paper, we shall con- 

centrate on simple ud-quark condensates and explore the 
possibility that these condensates are in one of the gapless 
phases with broken spatial symmetry [l3-[2l[. The ho- 
mogeneous gapless color superconductors are not treated 
because they develop tachionic instabilities in the ground 
state [2^ and their specific heat has anomalous temper- 
ature dependence [23j . 

Thermal components in the photon spectra from a 
handful of pulsars have been identified. The tempera- 
tures inferred from these measurements point toward a 
dichotomy in the evolution of isolated neutron stars - 
some of them are clearly much hotter than their same- 
age counterparts. The scatter in the data over the 
temperature-age plane could be correlated with certain 
static (mass, radius) and/or dynamic (spin, i?-field) pa- 
rameter(s) of pulsars. Such correlations are not easy to 
resolve theoretically, as there could be a number of fac- 
tors influencing the emitted spectrum (e.g., the prop- 
erties of the emitting regions could be sensitive to the 
variation of the magnetic fields at the surface of the star 
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and its composition). 

One of the aims of this work is to understand the cor- 
relation between the mass and the surface temperature of 
a cooling compact object for a class of models which sup- 
port high-mass {M ~ 2Mq) stars. Specifically, we shall 
study below the models constructed in Ref. [23] ■ The se- 
quences of static and rotating configurations of Ref. [23 ] 
have maximal masses of the order 2Mq with radii of the 
order of 12 km. The internal structure of these models is 
described in Ref. [25|, where it is shown that their mas- 
sive members contain quark cores with radii of the order 
of 7 km and masses < IM©. 

Previous work studied the cooling of quark-hadronic 
stars with different input physics and varying degree of 
realism. From simulations of cooling of compact stars 
with a mixed phase of quarks and hadrons, in which the 
neutrino emission processes are uniformly suppressed by 
a gap in quark spectrum, one may conclude that quark 
matter is invisible if the gap in quark matter is large 
and indistinguishable from hadronic matter if the gap 
is small p6j|. In simulations of the cooling compact ob- 
jects of 1.4 Mq with quark cores in the color- fiavor-locked 
(CFL) phase, the cooling is too fast to accommodate 
the observational data \27\. Because all the fermions in 
the CFL phase are gapped, the quark core remains hot 
for longer periods of time (up to 100 yr) compared to 
the phases containing ungapped fermions in the quark 
phase. The cooling simulations of Ref. of 1.4 Mq 
toy models having core quark matter in the gapless ver- 
sion of the CFL phase show that the cores of such stars 
would be much hotter than their nuclear counterparts 
in the photon cooling era. This is because such matter 
has larger specific heat due to quadratic rather than lin- 
ear dependence of the quasiparticle spectrum on the mo- 
mentum. The thermal behavior of Larkin-Ovchinnikov- 
Fulde-Ferrell quark matter was studied for 1.4 Mq toy 
models, and it was concluded that such phases will lead to 
rapid cooling of compact stars with rates that are conipa- 
rable to the rate found for the unpaired quark matter [29| . 
The quark stars that contain two-flavor quark matter in 
a homogeneous state, possibly in a mixed state with nor- 
mal quark matter, would cool too fast to account for the 
data [131. The possibility that the normal quarks addi- 
tionally pair in some channel fits the data [30, [3l| if the 
masses of the models are varied in a range 1.0 — 1.75 Mq 
and some suitable density dependence is assigned to the 
smaller gap. 

In this work, we are concerned with the thermal evo- 
lution of a sequence of compact objects, parametrized by 
their central density or mass, which contains large-mass 
(M < 2Mq) members with quark cores. Our models 
have realistic backgrounds, being derived from a micro- 
scopically motivated equation of state and the general 
relativistic structure equations. One aim of this study is 
to quantify the changes in the cooling behavior of stars 
along the sequence, focusing mainly on the new features 
that arise with the onset of quark matter in the cores 
of high-mass models. We assume quark matter of light 



u and d quarks in beta equilibrium with electrons; the 
s and other quark flavors are suppressed by their large 
masses. The pairing among the u and d quarks occurs in 
two channels: the red-green quarks are paired in a con- 
densate with gaps of the order of the electron chemical 
potential; the blue quarks are paired with (smaller) gaps 
of order of keV, which is comparable to core temperature 
during the neutrino-cooling epoch [32n34j . The neutrino 
emission in the case of single-flavor^ but multicolor pair- 
ing is as in normal quark matter |35| and differs from the 
single-color case studied here. 

A homogenous red-green condensate, when in the gap- 
less BCS regime, develops tachionic instabilities, which 
manifest themselves in negative Meissner masses [2^. 
Furthermore, the thermodynamic quantities of a homoge- 
neous gapless superconductor are generally anomalous at 
low temperatures; in particular, for large enough asym- 
metries the gap disappears at a lower critical temperature 
and the specific heat experiences a jump [l^il^. Recent 
studies show that a realization of inhomogeneous super- 
conductors, the Fulde-Ferrell phase, is preferred in color 
superconducting quark matter to the homogenous BCS 
phase at sufficiently low temperature and weak coupling 
for arbitrary fiavor asymmetries f20| as well as when color 
and charge neutrality are enforced (2l| . Furthermore, 
the anomalies in the thermodynamic quantities are cured 
when the inhomogeneous phases are treated (this conjec- 
ture is rigoro usly proven for the so-called Fulde-Ferrell 
phase 



Motivated by these observations and to avoid anoma- 
lies related to the homogeneous phases, we shall assume 
a generic inhomogeneous state and model the neutrino 
emission and specific heat in a phenomenological manner 
(see Secini. For the red-green condensate, we shall use 
a parameterization of neutrino emissivity in terms of the 
gaplessness parameter, as proposed in Ref. [13] . Thus, 
the second aim of this work is to understand the system- 
atics of the cooling tracks in the temperature versus age 
plane when the "gaplessness" parameter, i.e., the ratio 
of the pairing gap for red-green quarks to the electron 
chemical potential is varied. The magnitude of the gap 
in the spectrum of blue quarks is another parameter that 
will be varied. We do not attempt to fit the data by fine- 
tuning the available parameters (the mass of the star, 
the gaplessness parameter, and the gap for blue quarks); 
rather, we would like to find general trends that could 
provide orientation once the observational data demands 
such fits. 

This article is organized as follows. In Sec. |TT1 we de- 
scribe the models of purely hadronic and quark-hadronic 
stars, the pairing patterns in quark matter, and the mi- 
croscopic input needed for studying their thermal evolu- 
tion. Section IIIII studies the thermal evolution of these 
models. Here, we discuss the key approximations used to 
evolve the models in time and discuss the variations in 
cooling behavior as the central density of the models as 
well as the pairing parameters are changed. Our conclu- 
sions are collected in Sec. IIVI Some preliminary results 
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were reported in Refs. [37H3' 



II. QUARK CORES AND NEUTRINO 
EMISSION PROCESSES 

A. Equation of state, composition, and pairing 

We start with a brief description of the models under- 
lying our study. The low-density nuclear matter equa- 
tion of state is based on a relativistic model of nuclear 
matter in a homogeneous core and an inhomogeneous 
low-density phase constituting the crust of the star (see 
Ref. 24] for details). The high-density matter equa- 
tion of state is derived from the effective Nambu-Jona- 
Lasinio four-fermion interaction model which incorpo- 
rates the pairing interactions among quarks. The two 
equations were matched by a Maxwell construction, i.e., 
the requirement that the pressures of the two phases at 
the point of transition to the deconfined quark phase 
should match at a fixed baryo-chemical potential. The 
high-density quark phase pairs the up (u) and down id) 
quarks. We will assume that strange quarks (because 
of their large mass) are not important at the densities 
relevant to our models. With this equation of state as 
an input, the models of low-density purely hadronic and 
high-density quark-hadronic stars were constructed by 
solving the general relativistic structure equations. The 
masses and radii of various components are shown for a 
set of models in Table ID It is assumed that the transition 
from the hadronic core to the crust occurs at O.Spo, where 
Pq — 2.8 X 10^^ g cm~'^ is the nuclear saturation density. 
In the isothermal-interior approximation described be- 
low, the temperature gradients are concentrated in the 
envelope of the star; the transition from isothermal inte- 
rior to the nonisothermal envelope occurs at the density 
Pm = 10^" g cm^"^ and the radius i?cr- Table |T] lists the 
values of the inner radius Rcr of the models. 

Figure [T] shows the density profiles of four models 
selected from the set above, whose thermal evolution 
will be studied in Sec. IIIIl One of these is a purely 
hadronic model with central density Pc,i4 — 5.1, where 
Pc,i4 = Pc/(10^'* g cm~'^). The remaining three models 
contain a quark matter core and have central densities 
Pc,i4 = 10.8, 11.8, 21.0. The model with p^^u = 21.0 has 
the maximum mass for the sequence of the stars defined 
by the equation of state above, i.e., for larger central 
densities the stable branch of equilibrium, nonrotating, 
configurations terminates. The quark-hadron transition 
is seen in the density jump at the radius at which the 
phase transition occurs. (Note that the pressure and the 
baryonic chemical potential are continuous at the point 
of transition by construction.) The chosen models pro- 
vide a good coverage of the possibilities provided by the 
underlying equation of state. In particular, a substantial 
fraction of the mass and volume of the model with largest 
central density is occupied by quark matter. 

The combination of charge neutrality and beta equi- 



TABLE I: Masses and the radii of models constructed from 
the equation of state adopted in this work. The low-density 
models are purely hadronic, whereas the high-density models 
contain quark cores. The central densities of the models are 
given in units of lO'^* g cm"'', the radii in kilometers, and the 
masses in units of the solar mass. The columns from left to 
right are: the central density, the radius of the quark core, 
the quark-plus-hadronic core, the isothermal interior and the 
star, the masses of the quark core, the hadronic core, the 
crust, and the star. 



Pc,14 Rq Rq + H Rcr R Mq Mh Mc: 



M 



3.6 - 9.48 13.20 13.80 

5.1 - 11.47 13.39 13.53 

8.2 - 12.57 13.55 13.57 



0.43 0.11 0.54 
1.03 0.07 1.10 
1.81 0.04 1.85 



10.8 0.68 12.54 13.49 13.5 0.001 1.81 0.040 1.91 

11.8 3.41 12.40 13.31 13.32 0.09 1.80 0.04 1.93 

15.0 5.84 11.86 12.54 12.55 0.52 1.46 0.03 2.01 

21.0 6.77 11.34 11.91 11.92 0.89 1.14 0.02 2.05 



librium implies the presence in the stellar matter of elec- 
trons, whose chemical potential obeys the condition p.^ = 
Pd — Pu - Because of the shift in the chemical potentials of 
the quarks of different flavors, the pairing pattern differs 
from the conventional Bardeen-Cooper-Schrieffer (BCS) 
superconductors by the fact that the paired fermions are 
drawn from different Fermi surfaces. Under these condi- 
tions the pairing is "gapless." i.e., there exist segments 
on the Fermi surface where the gap separating the ex- 
cited states of the system from the ground state disap- 
pears. The gapless phases could be either homogeneous 
(in the sense that the spatial 0(3) symmetry is intact) 
or inhomogeneous, in which case the symmetry is broken 
down to a certain symmetry subgroup. 



B. Neutrino emissivities and other input 

Quark matter consisting of two light flavors cools via 
the beta decay (Urea) reactions 

d^u + e-\-D, u + e^d + u, (1) 

where v and u stand for neutrino and antineutrino. For 
unpaired quarks and to leading order in the strong cou- 
pling constant as the emissivity of the process per quark 
color is given by po[ 
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FIG. 1: Dependence of the density of the stellar matter on 
the radial coordinate for static, spherically symmetrical stars. 
The central densities and radii of the stars can be read off 
from the respective axes. The observed jumps in the density 
profiles of high-mass stars are a consequence of the onset of 
the quark matter in the core. 



where G — G cos 9 is the weak couphng constant, 6 is the 
Cabibbo angle, pd, Pu, and Pe are the Fermi momenta of 
down quarks, up quarks, and electrons. Non-Fermi-liquid 
corrections to the Urea rate above [4l| were neglected in 
our study. 

The quark pairing modifies the temperature depen- 
dence of the process In the BCS-type superconduc- 
tors, the process is suppressed nearly linearly for T 
and exponentially for T <^ Tc, where Tc is the crit- 
ical temperature [42]. In gapless superconductors the 
emergence of the new scale J/i — {fid ~ y^u)!"^ leads 
to two essentially different possibilities ;36], which are 
distinguished by the value the dimensionless parameter 
C = A/(S/z, where A is the gap in the spectrum in the 
limit ^/i = 0. When ^ > 1, the entire Fermi surface is 
gapped, i.e., when the thermal smearing of the Fermi 
surface ~ T is much smaller than the gap, no excitation 
can be created out of the Fermi sphere. In this respect, 
the C > 1 case resembles ordinary BCS superconductors. 
When C < 1, particles can be excited from the gapless 
regions of the Fermi sphere. We shall use below the in- 
terpolation formula proposed in Ref. [Ill, which covers 
these asymptotic limits, i.e., 

e7(C;T<T,) = 2/(C)e;3, 

/(C) = ,1^, ^ (3) 

exp[(C-l)fi-l] + l 

where the factor two comes from the two colors of quarks 
involved in the pairing. Since the pairing pattern breaks 
the S'C/(3) color symmetry, one of the quark colors (say, 
blue) remains unpaired. Same-color pairing leads to 



much smaller gaps in the spectrum of the quarks; these 
could be as small as tens of keV, and therefore irrele- 
vant in the neutrino-cooling era where T > 100 keV in 
the core. The pairing in this case is BCS-type, and we 
shall simply parameterize the modifications of the Urea 
process on blue quarks as in the case of baryonic matter, 
i.e., 

4(T<r,) = e^exp(^-^). (4) 

Note that = due to the different number of colors 
involved in these cases. 

As in ordinary BCS superconductors, the pair corre- 
lations suppress the specific heat of quarks (for T <C Tc 
exponentially); at T = Tc, the specific heat jumps (in- 
creases) if the phase transition to the superconducting 
phase is of second order. To model the modifications of 
the specific heat, we adopt the fitted formula for the ratio 
of the specific heat in the superconducting cs and normal 
Cat phases 

cs . , f (127r/7)(27T)-3/2e-V(7r) Q < r < 0.3, 
CN I -0.244 + 0.255t + 2.431t2 0.3 < t < 1, ^ 

where t = T/T^ is the temperature in units of the critical 
temperature and 7 — 1.781. Equation ^ is applied to 
neutron, proton, and blue-quark condensates. Note that 
Eq. ([5]) applies strictly only to isotropic BCS condensates. 
The anisotropy of neutron P-wave condensate somewhat 
modifies this ratio [i^ . 

For the inhomogeneous red-green condensate, one 
needs to rcscale the critical temperature which changes 
with 

r,(C)c.T,oyi-^<5(C), (6) 

where ^l = {^id + M«)/2, Ao = A(C = 0), T^o = Tc{C = 
0), and (5(C) = {rid - nu)/{nd + n„). The fully gapped 
and gapless regimes behave differently; in the presence of 
gapless modes, the specific heat has a linear dependence 
on the temperature. A phenomenological way to model 
this behavior is given by the relation 

4^(C;r<T,) = ficy^, (7) 

where is the specific heat of red-green unpaired 
quarks, taken as that for noninteracting quarks [45| and 
Cg^ is the specific heat of pair-correlated quarks. To com- 
pute the specific heat of the electrons, we assume they 
form a noninteracting, ultrarelativistic gas. 

The hadronic matter in our models is composed of 
neutrons (constituting ~ 95% of baryonic matter con- 
tent), of protons (with abundances of ~ 5% of bary- 
onic matter) and electrons (whose number is equal to 
that of protons by virtue of charge neutrality). In the 
unpaired state, the dominant cooling agents for such 



abundances are the modified Urea processes and neu- 
tral current bremsstrahlung processes, included in stan- 
dard fashion. We use emissivities derived for the free 
pion-exchange model of strong interaction '46*1, but re- 
duce them by a constant factor 5 to account for short- 
range repulsive component of the nuclear force. (Note 
that medium modifications of pion dispersion can fur- 
ther modify these rates and the resulting cooling substan- 
tially; see Ref. [13]. However, no clear evidence of such 
modification were found in experiments so far). The cor- 
responding rates below are exponentially suppressed, 
e.g., for the modified Urea process we have 



e^modCT" < Tc) = e^modexp 



A„(T) + Ap(T) 



T 



(8) 



In addition, we include the process of electron 
bremsstrahlung on nuclei in the crusts (48j . Following 
the latter reference, we assume that the ions form a fluid 
throughout the crust (or the number of impurities is high) 
so that the emissivity scales as T^. If at some tempera- 
tures ions form a crystalline phase, the emissivities scale 
roughly as T*" (classical crystal) or as (quantum crys- 
tal) and the bremsstrahlung emissivity is parametrically 
suppressed [49| . 

In the superfluid/superconducting hadronic phases, 
the processes of pair-breaking can contribute to the 
cooling [s^l- The neutral vector current processes 
are strongly suppressed by multiloop processes (sH - ls^ j. 
whereas the axial-vector emission can be taken at one- 
loop level to a good accuracy [H^l- Therefore, only the 
axial- vector neutrino emission through the pair-breaking 
processes were included in the simulations. The density- 
dependent zero-temperature pairing gaps of neutrons and 
protons, that have been used in the simulations, are 
shown in Fig. [2] Their non-zero-temperature values can 
be obtained in terms of the zero-temperature values A(0) 
using the formula [isj 



A(r) 
A(0) 



1 



^-7r/(7r) 



< T < 0.5, 
V3.016(l-r)- 2.4(1 -t)2 0.5 < t < 1, 



(9) 



which reproduces the exact BCS result with an error of 
order of a percent. 



III. THERMAL EVOLUTION 

Our models were evolved in time, with the input de- 
scribed above, to obtain the temperature evolution of the 
isothermal-interior. The isothermal-interior approxima- 
tion is valid for timescales i > 100 yr, which are required 
to dissolve temperature gradients by thermal conduction. 
Unless the initial temperature of the core is chosen too 
low, the cooling tracks exit the nonisothermal phase and 
settle at a temperature predicted by the balance of the 
dominant neutrino emission and the specific heat of the 
core at the exit temperature. The thermal conduction is 




FIG. 2: Dependence of the neutron and proton pairing gaps 
on their Fermi momenta, (solid line - neutron gap [62| |. 
dashed line - F2 neutron gap [63l] , and dashed-dotted 

line -^5*0 proton pairing gap [641]). 



effective enough to erase thermal gradients at densities 
larger than 10^'' g cm~'^. The low-density envelope main- 
tains substantial temperature gradients throughout the 
entire evolution; the temperature drops by about 2 or- 
ders of magnitude within this envelope. The isothermal- 
interior approximation relies further on the fact that the 
details of the temperature gradients within the envelope 
are unimportant if we are interested only in the surface 
temperature. Models of the envelopes predict the scal- 
ing T* = gsh{T), where gs is the surface gravity, and 
h is some function which depends on T, the opacity of 
crustal material, and its equation of state. We shall use 
the fitted formula Tg = 1.288(T4./5,i4)0-455 (Note 
that this relation does not take into account either the 
potential strong magnetization of neutron star envelopes 
or the presence of light elements in their surface layers. 
The modifications due to these factors are discussed in 
Refs. [5i,[53]). 

In the isothermal-interior approximation, the parabolic 
differential equation for the temperature reduces to an 
ordinary differential equation. 



-L,{T)-L.,{Ts) + H{T), 



(10) 



where L^, and are the neutrino and photon luminosi- 
ties, Cv is the specific heat of the core, and the heating 
processes, which could be important in the photon cool- 
ing era, are neglected, i.e., H{T) — (see Refs. [ssl. [59| 
for a summary of these processes). The photon lumi- 
nosity is given simply by the Stephan-Boltzmann law 
= AnaR^T^, where a is the Stephan-Boltzmann con- 
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FIG. 3: Time evolution of the surface temperature of four models with central densities (in units of P14) 5.1 (solid line), 10.8 
(long-dashed line), 11.8 (short-dashed line), 21.0 (dotted line). The integral parameters of these models can be found in Table|Tl 
The observational data points in order of increasing age correspond to RX J0822-4247, IE 1207.4-5209, RX J0002+6246, PSR 
1706-44, PSR 0833-45 (Vela), PSR 0538-^2817 (dots, red online), PSR 0656-^14, PSR 1055-52, PSR 0633-^1748 (Geminga), 
RX J1856. 5-3754, RX J0720. 4-3125 (squares, blue online). The dots correspond to fits assuming H-atmospheres, squares - 
to black-body fits. The data is taken from Ref. [6511. The triangles show in order of increasing age the current upper limits 
on the surface temperatures of PSR J0205+6449 [6l| PSR J112439. 1-591620 [H, supernova remnants (SNRs) G093. 3-1-6.9, 
G315.4-2.3, G084.2-F0.8, and G127.H-0.5 68], and RX J0007.0-h73 The pulsar in Gas A is shown by the plus sign 
The upper two panels correspond to cooling when the red-green condensate has ^ — 0.9, i.e., is not fully gapped; the lower 
panels correspond to = 1.1, i.e., the red-green condensate is fully gapped. The left two panels correspond to evolution with 
negligible blue-quark pairing (A;, — 0); the right two panels show the evolution for large blue pairing Aj, — 0.1 MeV. 



stant. The neutrino luminosities L^{T) and the specific 
heat Cy in Eq. pUj) are spatial integrals of emissivities 
and local specific heats of constituents over the volume 
of the isothermal-interior. In nonsuperfluid matter, the 
temperature dependence of the luminosities and the spe- 
cific heat in Eq. ([TU]) can be factored out and the re- 
maining coefficients need to be computed only once for 
a given stellar model. If matter is superfluid, the tem- 
perature and density dependence of the pairing gaps (or 
critical temperatures) cannot be factorized. Hence, we 
carry out the spatial integrations over the volume of the 
isothermal interior at each time step to obtain the rele- 
vant coefficients, vifhich are now temperature dependent. 
Once the coefficients are known, Eq. (fTO|) can be ad- 
vanced in time, e.g., through the Runge-Kutta algorithm. 
The results of integration of Eq. ([TU]) are shown in Fig. [31 
where we display the dependence of the (redshifted) sur- 



face temperature on time. Each panel of Fig. [3] contains 
cooling tracks for the same set of four models introduced 
in Sec. IIIII the cooling tracks for the purely hadronic 
model (solid lines) are the same in all four panels. The 
panels differ in the values of microphysics parameters, 
which characterize the pairing pattern in quark matter. 
Specifically, the two panels in the left column correspond 
to the case where the blue-quark pairing is negligible [i.e., 
the pairing is on a scale much smaller than the smallest 
energy scale involved, typically the core temperature). 
The two panels in the right column correspond to the 
case where the gap for blue quarks is large, A;, = 0.1 
MeV. The panels in the upper and lower rows are dis- 
tinguished by the value of the C parameter. [We use the 
values C = 0.9 (upper row) and Q = 1.1 (lower row)]. 

Figure [J] shows the dependence of the neutrino lumi- 
nosities of contributing processes on the core temperature 
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T[K] 



FIG. 4: Dependence of total luminosities of contributing processes on the core temperature: solid lines - Urea and 
bremsstrahlung processes in hadronic phase(s), long-dashed line - electron bremsstrahlung in the crusts, short dashed and 
dashed-dotted lines - Urea process with red-green and blue quarks, respectively, dashed-double-dotted lines - pair-breaking 
processes in the hadronic phases, dotted lines - photon emission from the surface. The underlying model corresponds to the 
1.93 solar-mass model in Table H] 



for the star model with central density pc.i4 — 11.8. It is 
complementary to Fig. [3] in the sense that it allows us to 
identify the dominant cooling agent for each temperature 
and, hence, via Fig. [3l for each epoch. Clearly, the cool- 
ing mechanism with largest luminosity determines the 
rate of decrease of the temperature. Since no heating 
processes are included in our study, the asymptotic rate 
of the cooling is determined by the photon luminosity. 
The dominant cooling agent in the neutrino emission era 
depends on the microphysics input as explained in cap- 
tion of Fig m 

A general trend seen in each panel of Fig. [3] is that 
the heavier the star, the faster is its cooling by neutri- 
nos. As a consequence, the heavy stars enter the photon 
cooling stage later than the lighter ones. In the late pho- 
ton cooling era, the arrangement of the cooling curves 
is inverted, i.e., the heavier stars are "hotter" than the 
lighter stars. (The photon luminosity of a star depends 
on its radius. Therefore, the heavier stars, which are 
more compact than the light ones, have smaller photon 
luminosities than the light stars of the same surface tem- 
perature; however, the cooling dynamics dominates over 
this static effect.) The fast cooling of massive stars in 



the neutrino emission era is driven by the quark Urea 
process (see Fig. 2]). Since the neutrino luminosity in 
Eq. ([TU| is an integral quantity, it scales with the size 
of the quark core Ly cx [i?']^, and therefore increases 
with the size of the quark core. The scatter of cooling 
curves in the T — t plane suggests that the experimen- 
tal data could be naturally accommodated in the models 
that have quark cores whose size changes as the central 
density of the model (and the net mass) are changed. 
In such an approach, the dichotomy in cooling behav- 
ior of neutron stars is attributed to the different masses 
of these objects: the lighter ones not containing quark 
matter cool more slowly than the heavier stars having 
quark cores. Note, however, that other factors, e.g., sur- 
face magnetic fields or the envelope structure, could be 
responsible for differences observed in the surface tem- 
peratures of isolated neutron stars as well. 

Superimposed on this general trend are the variations 
that arise from different functional dependence of the 
Urea emissivity of quarks on the temperature under dif- 
ferent assumptions about the pairing pattern of quarks. 
We explore the variation of the parameters Ah and (. The 
fastest neutrino cooling occurs for A;, = and C = 0.9 



8 



lU 






/- 

d 




















- '■' ' /' / 




J 10-" 




/' ■/■■ 


: 


_ 




- 


lo'^ 


/ 






lo'" 









10* lo' 

T[K] 



FIG. 5: Dependence of total luminosities on the core temper- 
ature for a 1.1 Mq hadronic star. Solid line - the luminosity of 
two-nucleon processes (modified Urea and bremsstrahlung), 
dashed-dotted line - the bremsstrahlung luminosity in the 
crust, dotted line - the photon luminosity of the surface, cir- 
cles - neutron pair-breaking process in the crust, squares - 
neutron and proton pair-breaking processes in the core. 



(the upper left panels in Figs. [3] and |4]). Indeed, in this 
case the Urea process is unsuppressed for all colors of 
quarks; the red-green quarks radiate neutrinos at rates 
comparable to unpaired quark matter, since there exist 
gapless excitations near the Fermi surfaces for ^ < 1. The 
next-fastest cooling scenario corresponds to the choice 
Ab = 0.1 MeV and C = 0.9 (the upper right panel in 
Figs, inland Fig. S]); here, the contribution of blue quarks 
is suppressed for T < 10^ K, whereas the radiation from 
red-green quarks is unsuppressed as before. Therefore, 
the neutrino luminosity is reduced by about 1/3 from 
the case where all quark colors radiate at the rate of un- 
paired quark matter. The effect of a large blue-quark 
gap is clearly seen in Fig.|4]where the initially power-law 
drop of the luminosity with temperature changes to expo- 
nential. Let us turn to the case where red-green quarks 
are gapped over the entire Fermi surface, i.e., C > 1. 
Because their gap is substantially larger than the blue- 
quark gap, their neutrino luminosity is suppressed at a 
much larger temperature. (This corresponds to the case 
C = 1.1, which is illustrated in the two panels in the lower 
row in Figs. [3] and HI) In the case where = (lower 
left panel) , we find the third- fastest neutrino-cooling sce- 
nario, in which the neutrino emission is dominated by the 
blue quarks and is about 30% of the rate that would be 
emitted by quarks of all colors. Finally, the slowest cool- 
ing scenario emerges when A;, = 0.1 MeV and C = 1-1 
(lower right panel) . The neutrino emission of quarks of all 
colors is shut off below the critical temperature of blue 
quarks. The photon cooling from the surface and the 
bremsstrahlung in the crust are dominant close to this 



temperature, but the photon luminosity takes over as 
the temperature drops slightly further (see Fig. |4l lower 
left panel). 

Figure |4] offers some further insight into the relative 
importance of various cooling processes. Specifically, the 
neutrino emission from the hadronic core with Lcore oc 
in the unpaired matter is an unimportant source of cool- 
ing except in the slow cooling model, because of the expo- 
nential suppression due to the gaps in neutron and pro- 
ton quasiparticle spectra. The bremsstrahlung process 
in the crusts with Lcrust oc is numerically unimpor- 
tant except for a limited range of temperatures prior to 
the transition to the photon cooling era. The fact that 
over long periods of time a single cooling process is domi- 
nant allows one to deduce the slope of the cooling curves, 
provided the specific heat capacity of the fermions scales 
linearly with temperature cy — aT and the nontrivial 
temperature dependence arising from pairing effects can 
be neglected. Writing the luminosity of a process as 
L — bT^^^'^, Eq. (fTO]) can be integrated to obtain 

, = ,. + £.n«)-(i-qf), (11) 

where is the initial time and To ~ T{t()) is the initial 
temperature. The expression in the braces is unity for 
T{t) < To, so one finds the scaling t - T"". Thus, the 
neutrino-cooling time via quark Urea process scales as 
T^*, while the cooling time via photons scales as T~^/^ 
for the surface-core temperature relation adopted above. 

In is further instructive to examine the relative impor- 
tance of different hadronic processes in a purely hadronic 
I.IM0 mass star. As seen in Fig. [5l in the temperature 
domain 7 x 10^ < T < Tc the neutrinos are produced 
predominantly by the pair-breaking processes in the core. 
At lower temperatures the crust electron bremsstrahlung 
and, to some extent, the pair-breaking bremsstrahlung 
in the crustal superfluid are the dominant processes. At 
about T ~ 10^ K the neutrino emission era ends and the 
photon cooling from the surface of the star becomes the 
dominant cooling process. 

A comparison with the measured temperatures of neu- 
tron stars suggests that slow coohng models (Af, = 0.1 
and ( = 1.1) describe young, cool objects. The tempera- 
tures of the fast cooling models are well below the current 
measurements and upper limits. Our 1.1 Mq hadronic 
model fails to account for the data on old and hot stars. 
There are several physical factors that may influence this 
outcome: (i) the rates of the neutrino bremsstrahlung of 
electrons in the crust may be overestimated due to the 
our assumption that the ions form a fluid pLOj : (ii) if 
the crust bremsstrahlung is subdominant due to crys- 
talline nature of the lattice, the pair-breaking process are 
the most important cooling agents in a wider tempera- 
ture range and these depend sensitively on the values of 
not-well-constrained pairing gaps in the baryonic mat- 
ter; (iii) variations in the magnetic fields and the surface 
composition of the neutron stars. 
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IV. CONCLUSIONS 

Motivated partly by the recent observation of a mas- 
sive compact object in a binary system Q and partly by 
the recent development of sequences of stable massive hy- 
brid stars with realistic input equations of state [l^, [l^] , 
we have modeled the thermal evolution of compact stars 
containing quark cores. These sequences of stable stars 
allow for a transition from hadronic to quark matter in 
massive stars {M > l.SSM©) with the maximal mass of 
the sequence ~ 2Mq. 

Our study has concentrated on the changes in the cool- 
ing behavior of the compact objects as the central density 
(and therefore macro-parameters such as the mass) is var- 
ied in a wide range along a sequence which has low-mass 
(M < 1.85Mq) purely hadronic members and high-mass 
{M > I.85M0) stars containing quark cores. We have 
seen that different cooling scenarios arise depending on 
the assumption of the magnitude of the pairing gaps in 
quark matter, including the possibility of emergence of a 
gapless phase in u-d condensate of quarks. 

We find that (i) the neutrino-cooling is slow for 
hadronic stars and becomes increasingly fast with an in- 
crease of the size of the quark core, in those scenarios 
where there are unpaired quarks or gapless excitations 
in the superconducting quark phase. The temperature 
scatter of the cooling curves in the neutrino cooling era 
is significant and can explain the observed variations in 
the surface temperature data of same age neutron stars, 
(ii) If quarks of all colors have gapped Fermi surfaces, 
the neutrino cooling shuts off early, below the pairing 
temperature of blue quarks; in this case, the tempera- 
ture spread of the cooling curves is not as significant as 
in the fast cooling scenarios, (iii) As the stars evolve 
into the photon cooling stage the temperature distribu- 
tion is inverted, i.e., those stars that were cooler in the 
neutrino-cooling era are hotter during the photon cooling 
stage. 

Before closing, we briefly state some potential modifi- 
cations to the results above and some further refinements 
of the models presented above. In this work, we concen- 
trated on the phase with two-flavor pairing. If the in- 
medium masses of strange quarks are low, they could ap- 
pear at low densities relevant to compact stars; the cool- 
ing of objects containing condensates with strange quarks 



has been studied elsewhere [28|. Since the fermionic ex- 
citations in these denser phases differ from those in our 
model, their cooling significantly differs from the cooling 
behavior we observe. 

The late-time evolution of compact objects will be af- 
fected by the heating in their interiors, as indicated in 
Eq. (ITO)) . Heating processes change the thermal history in 
the photon cooling, implying a less steeper than cx T~^/^ 
slope of temperature decay with time (ssl . [soj . While im- 
portant, these processes are unlikely to affect the present 
results and conclusion, which are mainly concerned with 
the neutrino emission era. 

Although our models have realistic density profiles de- 
rived from microscopic equations of state, we have as- 
sumed constant pairing gaps throughout the color su- 
perconducting phases. This is not the case, since the 
density of states at the Fermi surface and running cou- 
pling of strong interactions (both of which control the 
strength of pairing interactions) vary with density. As 
the density increases, the density of states increases as 
1/3 power of density, whereas the strong coupling con- 
stant decreases logarithmically. These variations would 
be of minor importance unless the C function occurs to 
cross C ^ 1 boundary separating the gapless and gapped 
phases. If this occurs, the core shell where ( < 1 will be 
the predominant neutrino emitter; physically, instead of 
the entire core volume the effective radiation volume will 
amount to that of the shell where C < 1- For any given 
model, this will reduce the neutrino luminosity of red- 
green quarks. A density-dependent ^ parameter, com- 
bined with other parameters discussed above, can be used 
to fine-tune the models to fit the observational data. In 
particular, our models can be applied to describe the re- 
cently observed cooling behavior of the neutron star in 
Cassiopeia A [iO, HH . 
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